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Abstract 

We discuss the mean-field approximation for a trapped weakly-interacting 
Bose-Einstein condensate (BEC) and its connection with the exact many- 
body problem by deriving the Gross-Pitaevskii action of the condensate. The 
mechanics of the BEC in a harmonic potential is studied by using a variational 
approach with time-dependent Gaussian trial wave- functions. In particular, 
we analyze the static configurations, the stability and the collective oscilla- 
tions for both ground-state and vortices. 

32.80.Pj; 67.65.+z; 51.10.+y 



1 



Typeset using REVIgX 



I. INTRODUCTION 



Few years ago, the Bose-Einstein condensation has been observed with alkali-metal va- 
pors confined in magnetic traps at very low temperatures. 1-3 By now, the total number of 
groups that have achived Bose-Einstein condensation is over a dozen. The confining poten- 
tial is accurately described by a harmonic trap, and theoretical studies of the Bose-Einstein 
condensate (BEC) in harmonic traps have been performed for the ground state, elementary 
excitations and vortices (for a review see Ref. 4). 

In the first part of this paper, we discuss the time-dependent Gross-Pitaevskii (GP) 
equation, 5 ' 6 which describes the macroscopic wave-function of a weakly-interacting Bose 
condensate at zero temperature. We deduce the GP action and the GP equation by imposing 
the least action principle to the quantum many-body action of the system within the mean- 
field approximation. 7 Such approach is mathematically equivalent to the field theory one, 
but it avoids the tricky issues related with spontaneous symmetry breaking. 

The rest of the paper is devoted to the study of BEC by using Gaussian trial wave- 
functions to minimize the GP action. Minimizing the GP action with a Gaussian trial 
function makes it possible to obtain approximate solutions for the condensate wave function, 
which can be treated, to a large extent, analytically. This method has been used by several 
authors to study the BEC static properties and the dynamics near the ground-state. 8-14 Here 
we extend its application also to vortex states. The presence of vortex states is a signature 
of the macroscopic phase coherence of the system. Moreover, vortices are important to 
characterize the superfluid properties of Bose systems. 7 

II. VARIATIONAL PRINCIPLE AND MEAN-FIELD APPROXIMATION 

Let us consider a iV-body quantum system with Hamiltonian H. The exact time- 
dependent Schrodinger equation can be obtained by imposing the quantum least action 
principle to the action 
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s = I dt < m\foQj - H\m > , a) 

where tp is the many-body wave-function of the system and 

< i/j(t)\A\TP(t) >= J d 3 r 1 ...d 3 r N TP*(r 1 ,...,r N ,t)Ai;(r 1 ,...,r N ,t) , (2) 

for any quantum operator A. Looking for stationary points of S with respect to variation 
of the conjugate wave-function ip* gives 

ih^ = , (3) 

which is the many-body time-dependent Schrodinger equation. 

As is well known, except for integrable systems, it is impossible to obtain the exact 
solution of the many-body Schrodinger equation and some approximation must be used. 

In the mean-field approximation the total wave-function is assumed to be composed of 
independent particles, i.e. it can be written as a product of single-particle wave-functions 
<f>j. In the case of identical fermions, ip must be antisymmetrized. 7 By looking for stationary 
action with respect to variation of a particular single-particle conjugate wave-function <p* 
one finds a time- dependent Hartree-Fock equation for each (pf 

where h is a one-body operator. In general, the one-body operator h is nonlinear. Thus the 
Hartree-Fock equations are non-linear integro-differential equations. Note that for alkali- 
metal vapors, which are quite diluite, the independent-particle approximation gives reliable 
results. 4 Obviously, it is possible to go beyond the independent-particle approximation by 
including correlations in the many-body wave-function. 8 This is particularly useful to repro- 
duce the properties of the superfluid 4 He, that is a strongly interacting system. 9,10 



III. GROSS-PITAEVSKII ACTION 



In this section we discuss the zero-temperature mean-field approximation for a sys- 
tem of trapped weakly-interacting bosons in the same quantum state, i.e. a Bose-Einstein 
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condensate. 7 In this case the Hartree-Fock equations reduce to only one equation, the Gross- 
Pitaevskii equation, 5 ' 6 which describes the dynamics of the condensate. As previously dis- 
cussed, this equation is intensively studied because of the recent experimental achievement of 
Bose-Einstein condensation for atomic gasses in magnetic traps at very low temperatures. 1-3 
The Hamiltonian operator of a system of N identical bosons of mass m is given by 

N *2 i N 

H = Y,(-^ + V ( ri )) + - £ V(r h r,) , (5) 

where Vo(r) is an external potential and V(r, r') is the interaction potential. In the mean- 
field approximation the totally symmetric many-particle wave-function of the Bose-Einstein 
condensate reads 

V>(ri, r N , t) = (j>(r u t) ... <j>(r N , t) , (6) 

where 0(r, t) is the single-particle wave-function. Note that such factorization of the total 
wave-function is exact in the case of a non-interacting condensate. The quantum action of 
the system is then simply given by 

d 



S GP = N Jdt< <f>(t)\ih- - h s \cf>(t) > , (7) 



where 



^ = -^V 2 + \/ o (r) + -(iV-l)|rfV|0(r')|V(r,r'). (8) 

We call Sqp the Gross-Pitaevskii (GP) action of the Bose condensate. By using the quantum 
least action principle we get the Euler-Lagrange equation 



ih-4>{v,t) = [-_V 2 + \/ (r) + (iV-l) / rf 3 rV(r,r')|0(r',t)| 2 ]0(r,t) , (9) 

which is an integro-differential nonlinear Schrodinger equation. Such equation and the effect 
of a finite-range interaction have been analyzed only by few authors (see, for example, Ref. 
14-16). In fact, at low energies, it is possible to substitute the true interaction with a 
pseudo-potential 7 



V(r,r') =5%-r') , (10) 

where B = A7ih 2 a s /m is the scattering amplitude and a s the s-wave scattering length. In 
this way one obtains the so-called time-dependent GP equation 5 ' 6 
f) fi 2 

^-0(r, t) = [- — V 2 + V (r) +B(N- 1) |0(r, t) | 2 ] 0(r, f ) , (11) 

that is the starting point of many calculations. 4 Note that the GP equation is accurate 
to describe the condensate of weakly-interacting bosons only near zero temperature, where 
thermal excitations can be neglected. At finite temperature, one must consider a generalized 
GP equation plus the Bogoliubov equations, which describe the quasi-particle elementary 
excitations of the condensate. 7 

IV. TRIAL WAVE-FUNCTION OF THE CONDENSATE 

In this section we use a Gaussian trial wave-function to describe the ground-state of the 
Bose condensate and to minimize the GP action. 

We consider a triaxially asymmetric harmonic trapping potential of the form 

V (r) = \mujl (uj\x 2 + u 2 2 y 2 + u 2 3 z 2 ) , (12) 

where u>i (i — 1, 2, 3) are adimensional constants proportional to the spring constants of the 
potential along the three axes. In the rest of this section we write the lengths in units of the 



characteristic length of the trap a = yh/(mcuo), the time in units Uq 1 , the action in units 
h and the energy in units Huj . By using the delta pseudo-potential to model the interaction 
between particles, the GP action is explicitly given by 

S GP = N J dt dh <P(r,t)[i^ + \v 2 - V (r) - l -B(N - l)|0(r)| 2 ]0(r,t) , (13) 

where B = 4na s /ao in our units. 

We want minimize this GP action Sqp by choosing an appropriate single-particle trial 
wave-function. A good choice is the following 



= ^mm^' i Mr , {-A +m ^} • (14) 

with (xi,X2,xs) = (x,y,z). Oi and $ are the time- dependent variational parameters. The 
Oi are the widths of the condensate in the three axial directions. Note that, in order to 
describe the time evolution of the variational function, the phase factor i/3i(t)xj is needed. 11,12 
The choice of a Gaussian shape for the condensate is well justified in the limit of weak 
interatomic coupling, because the exact ground-state of the linear Schrodinger equation with 
harmonic potential is a Gaussian. Moreover, for the description of the collective dynamics of 
Bose-Einstein condensates, it has already been shown that the variational technique based 
on Gaussian trial functions leads to reliable results even in the large condensate number 
limit. 15 - 17 

By inserting the trial wave-function, after spatial integration, the GP action becomes 

S GP = N J dtL^Pi,^) , (15) 

where the effective Lagrangian is given by 

L = ~[t M + + 2^ + + - JL-] , (16) 

with g = (2/7r) 1 / 2 (a s /ao) and g = g(N — 1). Again, by imposing the least action principle, 
we find six Euler-Lagrange equations 



A = ~TT~ » ' & i + ufa = ^3 + , (17) 



with % = 1,2, 3. Note that the time dependence of $ is fully determined by that of <7j. 
Moreover, for g = the six equations are exact solutions of the time- dependent GP equation. 
The last three differential equations are the classical equations of motion of a system with 
coordinates cxj and total energy per particle given by 

E 



= l[E\ d i]+ u fa>°*>°*)> (18) 

i=l 

where the potential energy is given by 



U(^a 2 ,a 3 ) = l[±(^ + ^) + -^-} . (19) 
2 l ~{ 2 lal en 02 03 J 



We use the previous equations to calculate the ground-state energy and the low-energy 
excitations of the condensate. 

The ground-state energy per particle of the Bose condensate is simply E^/N = U(a*), 
where a* = (er*, a^, cr%) is the the minimum of the effective potential energy U. This equi- 
librium point is identified by the following equations: 

24~°"l i 24- °2 , 2 4 ~ °3 / on \ 

- 9 = 1 > ^2 ff 2 - 9 = 1 , u 3 °3 - 9 = 1 • (20) 

a 2 a 3 (T1O-3 cricr 2 

The low-energy collective excitations of the condensate are the small oscillations of vari- 
ables aVs around the equilibrium point. The calculation of the normal mode frequencies for 
the motion of the condensate is reduced to an eigenvalue problem for the Hessian matrix A, 
given by 

d 2 U 



(21) 



1 /2 

In the case g = 0, the stationary solution is simply a* — (i — 1, 2, 3) and the energy 
per particle of the non-interacting condensate is given by E®/N = (wi + uo 2 + w 3 )/2. The 
Hessian matrix A of the potential energy U is diagonal and the frequency modes are simply 

= 2ui (i = 1,2,3). 

The general solution of Eq. ( |20"D for the stationary state and of the eigenvalue problem 
for the collective mode frequencies can be obtained numerically, 17 but analytical results are 
available in the Thomas-Fermi limit. In such limit, i.e. when g ^> 1, the right-hand side in 
Eq. (|20|) , related to the kinetic pressure, can be neglected and the following expression for 
the equilibrium point is obtained 

^ = \-^r) ' a * = \-^r) ' °* = y-^r) ■ (22) 

The energy per particle of the ground-state results 

— = 2(9^1^2^3) , (23) 
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and the Hessian matrix A is given by 



A 



/o 2 \ 

UJ\UJ'2 3 Uj\ 



(24) 



U2UJ3 3 u;J 

Note that 5 does not appear in the matrix. The eigenfrequencies Q of the collective motion 
are found as the solutions of the equation 14 

Q 6 - 3 (ul + w\ + oofj VL A + 8 (u\ w\ + cj 2 w 2 + w\ ufj Q 2 - 20 uol oo 2 2 uo 2 z = . (25) 

This equation has been also obtained by Stringari 4 using a hydrodynamical approach, which 
is effective in the Thomas- Fermi limit. 

For an axially symmetric trap, where uo\ = 002 = 00 j_ 7^ CJ3, the previous equations give 

a{ = a* = (^3/ci) 1 / 5 , 0$ = (Swj>f) 1/5 , £( )/7V = (5/4)^/ 5 4 /5 „ 3 2/5 , and 



fii, 2 = ^2cu 2 + p| ± ^16^1 + 9cu 4 - 16c^ 3 2 , fi 3 = cu ± . 



(26) 



Because of the axial symmetry of the trap, the angular momentum along the z axis is 
conserved and one can thus label the modes by angular (/) and azimuthal (m z ) quantum 
numbers. One finds 4 ' 11 m z = for the coupled monopole (/ = 0) and quadrupole (I = 2) 
modes f2i i2 , and \m z \ = 2 for the quadrupole (/ = 2) mode f2 3 . Note that the experimental 
results obtained on sodium vapors at MIT are in excellent agreement with these theoretical 
values. 18 

For an isotropic harmonic trap with frequency u, we obtain a* = {g/oo 2 ) 1 ^ (i = 1,2,3), 
E^/N = (5/4)(? 2 / 5 a// 5 , and ft lj2 = V5 u, Q 3 = V2 u. Here, differently from the axially 
symmetric case, the frequencies do not depend on m z . In particular, \fhuj is the monopole 
oscillation, also called breathing mode, characterized by radial quantum number n r — 1 and 
angular quantum number I = 0. Instead \/2u is the quadrupole (I = 2) surface (n r = 0) 
oscillation. 

It is worth noting that the variational formalism makes it possible to investigate the 
case of negative scattering length. In this case, the large N limit of the Thomas-Fermi 



approximation does not apply because the condensate number cannot exceed a critical value. 
Let us consider again the isotropic trap with frequency uj. For the ground-state, the number 
of particles satisfies the following equation 



where Cq — o* (i — 1, 2, 3) is mean radius of the condensate, that is the same in the three 
axial directions. It follows that, for g > 0, with N = 1 we have a = 1 and with N — > oo 
then do — > oo. Instead, for g < 0, the mean radius <To decreases by increasing the number 
of bosons iV to a critical radius erg = (5c<j 2 ) _1 / 4 with a critical number of bosons given by 
iV c = 1 + \g\~ l u~ l l 2 {h~ l l i — 5~ 5 / 4 ). Moreover, one finds a simple expression for the monopole 
(/ = 0) frequency of surface oscillation of the condensate, namely 16 



We have seen that a = 1 for iV = 1 or g = 0. Moreover, when N > 1, it is <7o > for g > 
and <7q < a < 1 for g < 0. Regarding the frequency f2, one easily verifies that fl — > v^cu for 
(j — > oo and Q — > for a ~" ^ °"o- Note that, for the triaxial trap, at the critical point the 
two highest frequencies diverge and the lowest one falls to zero. 17 

It is important to observe that it is possible to include in the Gaussian trial wave-function 
other three parameters describing the position of the center of mass of the condensate. Due 
to the harmonic confinement, the motion of the center of mass (dipole mode) is periodic and 
the frequencies of oscillation are equal to those of the harmonic trap. 11 ' 12 



Let us consider states having a vortex line along the z axis and all bosons flowing around 
it with quantized circulation. The observation of these states would be a signature of 
macroscopic phase coherence of trapped BEC. The time- dependent variational approach 
can be used to describe also these vortex states of the condensate in a axially symmetric 
trap 




(27) 



Q = 5u 2 — a. 



■41 1/2 



(28) 



V. TRIAL WAVE-FUNCTION OF VORTEX STATES 
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V (r) = Uiuir 2 + LU 2 3 z 2 ) , (29) 



2 

where r = y/x 2 + y 2 , 9 =aictg(y/x) and z are the cylindrical coordinates. 

The trial wave function of the vortex states of the condensate can be chosen as 



(30) 

where k is the vortex quantum number and cr±, a 3 , and /?3 are the time- dependent 
variational parameters. This trial wave-function describes exactly the vortex states of a 
non-interacting gas of bosons in the harmonic trap. 

By following the same procedure of the previous section, namely by inserting the trial 
wave-function in the GP action, and after spatial integration, one gets 

S gp = nJ dt L0 ± , /3 3 , ± , p 3 , <r ± , <j 3 ) , (31) 



where the effective Lagrangian is given by 



2 



1 _ o ~n 1 



(k + l)(P ± a 2 ± + ^ + 2a 2 ± (3 2 + -ooW) + 



+ (/V 3 2 + + + + , (32) 

2cr§ 2 a\a 3 



with g = (2/7r) 1 / 2 (a s /ao) and g(k) = g(N - l)(2k)\/ (2 2k k\ 2 ). The four Euler-Lagrange 
equations of the system are given by 

"-—it- k =~k' (33) 

(fc + l)* ± + (fc+lVi«T ± =^±^ + 4^-, *3+^V 3 = ^ + 4^- (34) 

The time dependence of /3± and /3 3 is fully determined by that of cr± and a 3 ; moreover, 
for g = (non-interacting condensed vortex) these four equations are exact. The last two 
differential equations correspond to the classical equations of motion of a system with energy 
per particle 
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i 



(k + l)&] 



-CTn 



+ U(a ± ,a 3 ) , 



(35) 



where the potential energy is given by 

U (a±, a 3 ) = -[(& + lplal + + 



(k + 1) 



+ 



2al 



a\a 3 



(36) 



cr ± ^ 3 

The equilibrium point (cr^,^), corresponding to the minimum of the effective potential 
energy, is given by the following equations: 



(k + l)u 2 ± ai-g(k)- = (k + l), 

03 



^303-^(^)^1 = 1 • 
(5 i 



(37) 



In the non-interacting case (g = 0) the stationary solution is simply a]_ = uj 1 / 2 an d 
l li 

o\ = uj 3 . The stationary energy of the non-interacting vortex state is exactly given by 



E, 



(0) 



N 



(k + 1)cj_l + ^uj 3 



(38) 



It is important to observe that the energy grows only linearly with the vortex quantum 
number k, in contrast to the homogenous case 7 where the vortex energy grows as k 2 . In 
addition, the two (m z = 0) frequency modes are Qi = 2u± and f2 2 = 2uj 3 and do not depend 
on the vortex quantum number k. 

In the Thomas- Fermi limit (g » 1), the analytic expression of the equilibrium point 
reads 



cr = 



g{k) uj 3 



nl/5 



_(k + l) 3 / 2 ul 
The stationary energy of the vortex state results 



~g{k) (k + l)uj 

<j4 



nl/5 



(39) 



4 0) 5 



(40) 



N 4 

From this formula it is easy to calculate the critical frequency Cl c at which a vortex can 
be produced. One has to compare the energy of a vortex state in frame rotating with 
angular frequency Cl, that is — ClL z , with the energy E° of the ground state (with no 
vortices). 7 Since the angular momentum per particle is hk, the critical frequency is given by 
hCl c = {Ef ] /N - EW/N)/k. By using our previous formula, we get 
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n ^ 1 4/5 2/5 



~g{kfl\k + lf' b -m 



2/5 



\2/5 



(41) 



In the presence of a vortex, for frequencies smaller than the critical frequency Q c , the system 
is not thermodynamically stable (global stability), but can be locally stable (metastability) 
against small deformations of the system. Rohksar 19 has shown that the vortex state is 
locally unstable for Q = 0. Actually, in a harmonic trap, there is local instability for 
< Cl < Cl s , where the frequency Cl s (Cl s < Q c ) can be estimated numerically. 20 On the 
other hand, we have recently shown that, in toroidal traps, vortices can be locally stable 
also in absence of external forced rotation. 21 

It is important to observe that, in our present treatment, the condensate has a unique 
vortex line, which is fixed along the z axis. Butts and Rokhsar 22 have numerically shown 
that increasing the rotational frequency well beyond C will favour the creation of more 
complex vortical configurations, associated with the occurrence of 2 or more vortex lines. 

To calculate the collective oscillations in the Thomas- Fermi limit, we follow the same 
procedure of the previous section. The Hessian matrix A of the potential energy U at the 
equilibrium point is given by 



A 



1 8 (k + 1) uj 2 ± 2 (k + l) 1 / 2 uj ± uj 3 ^ 



(42) 



^2 {k + lf 2 LJ ± LU 3 

Again the g{k) dependence drops out. Moreover, the mass matrix of the kinetic energy is 
not the identity matrix but instead 



M 



f 2 (k+ 1) 0^ 



V 



(43) 



/ 



I 

The eigenfrequencies Q of the collective motion are found as the solutions of the equation 
| A — VL 2 M\ = 0, which gives the m z = frequencies 



Q 1;2 = \\2ul + -cu 2 3 ± -y/16ui + 9cul - 1Qlu 2 ± ujI . 



(44) 



This is an interesting result: The m z = frequencies of a harmonically trapped condensate 
do not depend on the vortex quantum number k. Note that Zambelli and Stringari 23 have 
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shown, by using a sum-rule approach, that the m z ^ collective frequencies do instead 
depend on k. Nevertheless, all the modes are k sensitive because of the amplitudes, which 
are given by < r 2 > 1 / 2 = {k + 1) 1/2 (J_l and < z 2 > l / 2 = a 3 /V2. 

Finally, we study an isotropic trap with frequency u. We have a* ± = only for g — 0. 
In order to obtain some analytic results also in the case of negative scattering length, we 
consider the vortex trial wave-function with a unique variational parameter a, i.e. we put 
a± = o" 3 = a. For the ground-state, the number of particles satisfies the equation N = 
1 + g(k)^ 1 (u 2 a^ — (To)(l + 2/3k), where o"o is mean radius of the condensate in the stationary 
configuration and g(k) = g(2k)\/(2 2k k\ 2 ). It follows that for g < 0, the critical number 
of bosons for the collapse of the wave-function is given by iV c = 1 + \g(k)\~ 1 uj~ 1 ^ 2 (5~ 1 ^ — 
5~ 5 / 4 )(l + 2/3/c). Thus, compared with the condensed state without vortices, the critical 
number of trapped atoms in the vortex state increases. This variational result is in full 
agreement with the k — 1 numerical results of Dalfovo and Stringari. 24 

VI. CONCLUSIONS 

We have studied a Bose-Einstein condensate trapped in a harmonic potential by using a 
variational method based on the minimization of the Gross-Pitaevskii action with a Gaussian 
ansatz for the shape of the macroscopic wave function of the condensate. Our analysis 
has been applied to the whole range of coupling strengths between particles (attractive 
and repulsive interaction) to study both the ground-state and vortex states. In particular, 
for vortex states, we have determined the critical rotational frequency at which a single 
vortex becomes globally stable. Moreover, we have calculated the collective oscillations 
and the critical number of bosons for the collapse of the condensate in the case of attractive 
interaction. Nowadays it is experimentally possible to trap and condense a very large number 
of alkali atoms (more than 10 9 ) and to obtain a less diluite system. Thus, the study of the 
time-dependent variational approach without the independent-particle approximation begins 
to appear an interesting problem also for Bose gases in an external potential (the case of a 
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superfluid film containing a vortex has been discussed in Ref. 9 and 10). This will be one 
of our future projects. 
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